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Abstract 

We study the pair production of neutral Higgs bosons through gluon fusion 
at hadron colliders in the framework of the Minimal Supersymmetric Stan- 
dard Model. We present analytical expressions for the relevant amplitudes, 
including both quark and squark loop contributions, and allowing for mix- 
ing between the superpartners of left- and right-handed quarks. Squark loop 
contributions can increase the cross section for the production of two CP-even 
Higgs bosons by more than two orders of magnitude, if the relevant trilinear 
soft breaking parameter is large and the mass of the lighter squark eigenstate 
is not too far above its current lower bound. In the region of large tan/?, 
neutral Higgs boson pair production might even be observable in the 4b final 

state during the next run of the Tevatron collider. 
14.80.Cp, 14.80.Ly 
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I. INTRODUCTION 



The electroweak symmetry breaking sector is the last ingredient of the Standard Model 
(SM) that remains to be explored experimentally. If the masses of the W and Z bosons, 
as well as of the charged leptons and quarks, are due to the vacuum expectation value 
of an elementary scalar Higgs field, naturalness arguments indicate [JTJ] that Nature should 
become supersymmetric at an energy scale not far above the weak scale; otherwise quantum 
corrections would destabilize the hierarchy between the electroweak scale and the much 
larger Planck, GUT or string scale. Supersymmetry not only demands the introduction of 
superpartners for all know SM particles, but it also requires the existence of at least two 
Higgs doublets (and their superpartners). In its minimal version, the so called minimal 
supersymmetric standard model (MSSM) 0], there are three physical neutral Higgs bosons 
after electroweak symmetry breaking. In this work we assume that CP is conserved in the 
Higgs and squark sectors of the theory, and consequently the neutral Higgs bosons can be 
classified as two CP-even states h, H (with < inn) and one CP-odd state A. 

The production and detection of neutral Higgs bosons at e + e~ colliders is quite straight- 
forward ||. It can be shown || that at least one MSSM Higgs boson would have to be 
found at an e + e~ collider operating at a center-of-mass energy > 250 GeV. However, 
there are no immediate prospects for the construction of such a collider, and the energy of 
LEP may not suffice to find even one MSSM Higgs boson. On the other hand, the Teva- 
tron will soon begin its next collider run, at a slightly increased energy (^/i = 2 TeV) and 
greatly increased luminosity (anticipated / Cdt = 2 fb _1 per experiment); future runs with 
yet higher luminosity are being contemplated. In a few years experiments at the LHC will 
also start taking data at = 14 TeV and luminosity between 10 and 100 fb _1 /yr. It is 
therefore important to explore all channels that might yield information about the Higgs 
sector at hadron colliders. 

The largest contribution to the inclusive neutral Higgs production cross section at such 
colliders comes from gluon fusion, which can produce a single Higgs boson through quark 
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|5||| or squark || loops. The next-to-leading order (NLO) QCD corrections (involving 
gluons and light quarks) to these processes have already been computed [0,|J. On the 
other hand, NLO SUSY QCD corrections (involving gluinos) have not yet been calculated, 
however, they are expected to be smaller. One drawback of this production channel is that 
it can give a visible signal only if the Higgs boson decays into a mode with a small branching 
ratio. For example, in case of the light CP-even scalar h, the only promising decay mode 
in this production channel is /i — > 77, which has a branching ratio of the order of 10 -3 [f|. 
Moreover, if squark loop contributions are large, they tend to be destructive fli~0l , at least 
for gg —>■ h — > 77. As a result, gg —>■ <p production in the MSSM (0 = h, H, A) cannot give 
a viable signal at the Tevatron, and is often difficult to detect even at the LHC. 

One alternative is to search for the production of a single Higgs boson produced in 
association with some other particle(s) ||11|| . However, the main (hadronic) decay channels 
of MSSM Higgs bosons seem to be detectable in this way only in Wh (Zh) production at the 
Tevatron |]12[ if more than 10 fb _1 of data can be accumulated, as well as in some cases for h 



bosons produced in the cascade decays of heavy sparticles at the LHC |T3| . The detectability 



of hadronic Higgs decays in tt(f> production |I4J still awaits confirmation by an experimental 



study; no background estimates for the recently suggested |T5| tt(f) channels exists as yet. 
Finally, for very large values of tan/3 associate Higgs bb production might also be observable 
at the Tevatron in the 46 |1| and/or the bbr + T~ [[L7| mode. 

In this paper we instead study the production of two neutral Higgs bosons. It is hoped 
that the 46 final state, with invariant mass peaks in both bb pairs, will give a detectable signal 
at the Tevatron and/or the LHC at least in some regions of parameter space. Moreover, the 
Higgs-boson pair production can in principle be used to probe the Higgs self-couplings in 
order to reconstruct its potential [0. This process has first been discussed in refs. (19| for 
the SM, and in []20| the quark loop contribution in the MSSM has been studied. NLO QCD 
corrections to heavy quark loops have also recently been calculated JHj]; they are of similar 
relative size as the corresponding NLO corrections to single Higgs production. Here we 
extend these analyses by including the contribution from squark loops, allowing for general 



mixing of the superpartners of left- and right-handed quarks. We find that, unlike for single 
Higgs production ||10|| , squark loop contributions can increase the total cross section by more 
than two orders of magnitude. In some regions of parameter space beyond the reach of LEP, 
the 46 final state might even give a detectable MSSM Higgs signal at the next run of the 
Tevatron. 

The primary purpose of this paper is to present all analytical expressions required for a 
calculation of the squark loop contribution to the production of two neutral Higgs bosons, 
and to illustrate their potential importance with a few numerical examples. 

The remainder of this paper is organized as follows. In Sec. II we give analytical ex- 
pressions for all one-loop amplitudes of the form gg — > 4>i<j)j, where <pij are neutral Higgs 
bosons with definite CP quantum numbers. These expressions are completely general; a list 
of the relevant couplings as predicted by the MSSM is given in the Appendix. In Sec. Ill 
some numerical results are shown. We focus on scenarios with either a small or a large 
value of the parameter tan/5, which governs the size of Yukawa couplings in the MSSM. We 
find potentially very large squark loop contributions in both cases, but the experimental 
discovery of Higgs boson pairs at the Tevatron seems to be somewhat more promising at 
large tan/3. Finally, Sec. IV is devoted to a brief summary and some conclusions. 



II. CROSS SECTIONS AND MATRIX ELEMENTS 

We write the subprocess (parton-parton) differential cross section as, 

<*» _ |M| 2 

di 167TS 2 ' 1 ' 

where the hatted Mandelstam variables refer to the parton-parton system. The total cross 
section can be obtained from ([!]) by integrating over t and convoluting with the gluon 
densities in the two colliding hadrons, as usual. In general, the invariant amplitude can be 
written as 

|M|2 = 2 (j) (s) * 
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Here M 



(n) 



is the helicity amplitude for the production of the Higgs boson pair 



{4>i4>j) for the initial gluon helicities Ai(2) = ±. The sum runs over all Feynman diagrams 
(n) that contribute to a specific process. The factors refer to the average over the initial 
gluon helicities (7), the color factor [Tr(T a T b )] 2 = |5 aa = 2, and the average over the gluon 
colors (gj). The symmetry factor, Sf, equals \ for the production of two identical Higgs 
bosons, and is 1 otherwise. 

The Feynman diagrams contributing to the gg — > ifif , /iif , and AA processes are 
presented in Fig. [I], while the contributions to the processes gg — > hA and iJA are shown 
in Fig. |[ We now list the resulting helicity amplitudes for these two classes of processes, 
starting with the case where both produced Higgs boson have the same CP properties. We 
have used the FeynCalc package p2| for the analytical calculation. 



A. Invariant Amplitudes for gg — > hh, HH, hH, and AA 



1. Quark Loop Contributions 



By assumption these processes conserve CP, which implies that M ++ = M and M + _ = 

M_ + . The two independent helicity amplitudes for the production of two CP-even Higgs 
bosons HiHj = 1,2 for h and H respectively), where the superscript on the amplitudes 
refers to the number of the Feynman diagram in Fig. 1, are as follows. 



M^l(HiHj) = 



7T 



' "• [2+ (Am q 2 - s)C{0,0,s 



[(s-m Hk 2 )+im Hk T Hk 



(3a) 
(3b) 



MSW) = ^V Hm V Hm {-4S - 8m 2 q C(0, 0, s)s - (8 



™\ - m 2 H . - m 2 H . 



+ U i C(m% i Au)+T j C(m 2 H .,0,i) + U 3 C(m 2 H .0, 



(m H m H - tu)D(m H . , 0, m H , 0, t, u) 



2m 2 (8m 2 — m 2 H _ — m 2 H . — s)s x 



(3c) 



M {2) _{H i H j 



D(m 2 H .,0, m 2 H .,0, i, u) + D(m 2 H ., m 2 H ., 0, 0, s, i) + D{m 2 H ., vr? H . , 0, 0, s, &)] } 
ia^V^ { (8m 2_ f _^ )(2<m2ffj | _f._ fl2)C 7( m ^m^ l S) 



2ir(m 2 H .m 2 H . — tu) 



+ (m 2 Hi m 2 H .-8m 2 q t + t 2 ) x 



T.CV! . 0, t) + r i C(m| r ., 0, t) - §C(0, 0, 5) + stD{m 2 m 2 ,0, 0, 5, t) 



+ ( m Hi m % ~ 8m lu + U 2 ) X 



UiC(m 2 Hi ,0,u) + UjC(m 2 H .,0,u) -sC(0,0,s) + suD(m 2 H .,m 2 H ., 0, 0, s, u) 



+ 2mJ(m^mi. - tu)(8mj -t-u)x 



(3d) 



£>(m#., 0, m^, 0, i, u) + D(m 2 H .,m 2 H .,0, 0, s, i) + D(m 2 H .,m 2 H .,0, 0, s, u)] }3e.) 

Here, V# fcg9 and VHiHjH k are the Yukawa coupling constants of -fffc to quarks and the trilinear 
Higgs couplings, respectively. Expressions for these couplings in the MSSM are given in 
Appendix A while we list in Appendix B our choice for the polarization vectors. We have 
also defined the quantities 



Ti = (m 2 H . — i) and Ui = (m 2 H . — u) . 



(4) 



Furthermore, the loop integrals appearing in Eqs. (§) are defined in terms of the Passarino- 



Veltman scalar functions Co and Dq p3| (see Appendix B for our conventions) as 



C(a,b,c) = C (a,b,c,m g ,m q ,m q ) ; 
D(a, b, c, x, y, z) = D (a, b, c, x, y, z, m q , m q , m q , m q ) 



(5a) 
(5b) 



The scattering amplitudes for the production of a pair of pseudo-scalar Higgs boson are 



M { + ] _(AA) = 



(6a) 
(6b) 
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M^l(AA) = ^V Aqq V Aqq {2s + 4m 2 q C(0, 0, s)s 

ITS L 



- 2m l 



T A C(m A , 0, i) + U A C(m A , 0, u) + m A (m A — ui)D(m A , 0, m 2 A , 0, i, u) 



m 2 s(i + u) X 



(6c) 



D(m A , 0, m^, 0, t, m) + D(m A , m A , 0, 0, s, i) + D(m A , m A , 0, 0, s, u) } (6d) 
-V Aqq V Aqq {s{2m\ + f 2 + u 2 )C(0, 0, s) - 2T A (m A + P)C(rr4 0,t) 



— *a. 



27r(m 4 i — 



- 2U A (m\ + u 2 )C(m 2 A , 0, u) + (t + w) (2m^ - t 2 - u 2 )C{m% m% s) 
+2m 2 (i + u) (m 4 A - iu) x 



(6e) 



— s 



D(m A , 0, m A , 0, i, u) + D(m A , m A , 0, 0, s, t) + D(m A , m A , 0, 0, s, u) 
i{m A + i 2 )D(m 2 A} m A , 0, 0, s, t) + u(m\ + u 2 )D(m A , m A , 0, 0, s, «)] } (6f ) 



Analogously to the production of CP-even Higgs bosons, we have defined the quantities 



T A = (m A — t) and U & 



[m A — u ) 



(7) 



As a check of our calculations we verified that our results for the quark (squark) loop 
contributions to the Higgs pair production are invariant under QCD gauge transformations. 



Furthermore, our results agree with those of ref. |20 . 



2. Squark Loop Contributions 

We now turn to the new results of this paper, i.e. the squark loop contributions depicted 
in diagrams (3)-(8) in Fig. 1. These can be grouped into three sets of diagrams, (3) + (4), 
(5) + (6) and (7) + (8), which are gauge invariant and finite. Moreover, diagrams (7) and (8) 
are finite by themselves and we therefore list their contributions separately, treating QCD 
interactions in the 't Hooft-Feynman gauge. 

Vi 
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^sVHiqkqkVHiHiHj 
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l + 2miC kkk (0,0,s) 



(8a) 
(8b) 
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2vr 



-V, 



1 + 2m? C fcfcfc (0, 0, s) 



M^ 6) (H l H J ) = ; 



M++{HiH 3 ) 



IT 



M^iHiHj) = ; 



V 

2 n g H iQkm VH iikqi 



x lT\C lkk (m 2 H .,0, t) + UiCkuimjj.A u) + TjCku(m 2 H .,0, t) 
+ U j Ci kk (m 2 H 0,u) +2sC k i k (m 2 H .,m 2 H s) + 



m \ ~ m \d § ~ ( m H i m H j - *«)] D ikki(m 2 H .,0, m 2 H .,0, i, u) 
+ 2sm? \Di m (m 2 H 0, mi., 0, t, u) + D k i kk (m 2 H m 2 H 0, 0, s, t) 



+D klkk (m 2 Hi , m 2 H , 0, 0, s, u) 



MXf. 



2n(m 2 H .m 2 H - tu) 



^Hiq k qiVH.jq k qi 



x <( s(2m? fe - 2m? + t + u)C kkk (0, 0, s) 



9; 



-f 



U 



TiCi kk (m 2 Ht , 0, t) + TjCi kk (m 2 Hj ,0, i) 
UiCi kk (m 2 Hi ,0,u) + UjC lkk {m 2 H ,0,u) 



TAm 



n"«q k 



C kU (m 2 Hj ,0,t) + Ci kk (m 2 Hj , 0, t) 

- ^("4 - m |) [CfciK^HiJ °> «) + c ikk(m 2 Hi ,0, u) 
+ (2m|.m|. - t 2 - m 2 ) C Mk (m 2 Hi ,m 2 Hj ,s) 

+ ( m l - m l) 2 + ( m L + m l) - fe ) 

x [Afcfci( m fl-i> °> m Hj-> °> *> «) + D kikk(rn 2 H .,m 2 H .,0, 0, s, i) 
+D klkk (m 2 H . , m 2 ^ , 0, 0, s, u) 
+ [-s* 2 - (m| - m? ) (2ts - {m 2 H m 2 Hj - *«))] D klkk (m 2 H .,m 2 H .,0, 0, s, 
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-su 



III 



<lk 



D klkk (m 2 H m 2 H 0,0, s,u) 



MSSM predictions for the trilinear Higgs-squark-squark couplings VH t q k q k and the quartic 
Higgs-Higgs-squark-squark couplings V^H^qk are given in Appendix A; note that we need 
the former two couplings only for two identical squarks. The loop functions appearing in 
Eqs. (|]) depend on the squark masses and are defined as 

C ijk (a, b, c) = C (a, b, c, m Qi , m qj , m qk ) ; (9a) 
D ijkl (a, b, c, x, y, z) = D (a, b, c, x, y, z, m qi ,m q .,m qk ,m qi ) . (9b) 

The corresponding expressions for a pair of pseudo-scalar Higgs bosons are: 



M^\AA) = M^\H % HA\V HlH%Hj -> V HlAA \ ; 
M^ i] (AA) = ; 



r(3+4) 



M^ e \AA) = M^'^HjWh^ 
M^+ 6 \AA) = ; 
M^ + (AA) = M^iH.HjW^..^ -> V AMl ] ; 
M { +l{AA) = ; 

M^l(AA) = -M^H^Vh^ -> V Aqkqi ] ; 
M^l(AA) = -M^H.HjWh^ -> V AMl ] . 



r(5+6), 



V 



AAq k q k \ i 



(10a) 
(10b) 
(10c) 
(lOd) 
(lOe) 
(lOf) 

(10g) 
(lOh) 



B. Invariant Amplitudes for gg —>■ Ah, AH 

We now turn to the production of two Higgs bosons with different CP quantum numbers 
which only receive contributions from quark loops. Since we assumed CP-invariance in the 
Higgs and squark sectors, the squark contributions to the production of a CP-even and a 
CP-odd Higgs boson, i. e. the diagrams (4), (5), (6), and (7) of Fig. add to zero. Let us have 
a closer look at this. Note that A only couples to two different squark mass eigenstates; 
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this immediately eliminates the equivalent of diagrams (5) and (6) in Fig. 1. Moreover, 
^Aqiq2 = ~ ^a<Md an d consequently the two possible orientations of the loop momentum in 
diagrams (6) and (7) in Fig. 2 exactly cancel each other. Finally a same-flavor qq* pair 
in a color-singlet state coupling to a Z boson is in a CP-even state, while one CP-even 
and one CP-odd Higgs boson coupling to a Z boson are in a CP-odd state; this explains 
why diagrams (4) and (5) vanish. For completeness we give expressions for the quark loop 



contribution in our notation; our results agree with those of ref. |[20| . Note that we again 

only have two independent helicity amplitudes since M ++ = —M and M + „ = — M_ + for 

Ah and AH production. 



M^ + (AH- 



ia s m q 



7r[(s - m A 2 ) + im A T A ] 



V Aqq V„. AA s C(0, 0, s) 



M^l(AHi) = 



MfUAHi 



ta s g z T«(m 2 -m 2 A )(s-M 2 ) 



27im 2 z [(s — mz 2 ) + imzT z ] 



Vi 



ZAHi 



1 + 2m 2 C{0,0,s) 



M (2) _(AH l ) = ; 



V Hi qqV Aqq 



2-ks 

x <! (rn 2 A -m 2 Hi ) 



T A C(m 2 A , 0, t) + T % C{m 2 , 0, t) + U A C(m 2 A , 0, u) 



+ UiC(m H .i 0, u) — {rn A m H . — tu)D{m A , 0, rn H _, 0, i, u) 



+ 2m 2 s(T A + U A ) D(m\,0,mjj.,0,t,u) + D(m\,m%.,0,0,s,t) 



+D(m 2 A ,m 2 H ,0,0, s,u) 



M^_(AHi 



—la. 



-VH iqq V Aqq 



2ir(m 2 A rn 2 H . — tu) 

x -u 2 )C(0,0,s) + (Am 2 A m 2 Hi - (i + u) 2 )(i - u)C(m 2 A ,m 2 Hi , s) 

+ (rn 2 A m 2 H . - i 2 ) [T A C(m% 0, t) + T,C(m^, 0, t) 
+stD(m A , m 2 H ., 0, 0, s, i) 



(11a) 
(lib) 
(11c) 
(lid) 



(lie) 
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— (m 2 A m 2 Hi — u 2 ) U \C{m 2 A , 0, u) + UiC(m 2 H ., 0, u) 

+suD(m 2 A , rriff., 0, 0, s, u) 
+ 2m 2 q {t - u)(m 2 A m 2 Hi - iu) [D(m 2 A , 0, m 2 H ., 0, t, w) + D(m A , m 2 Hv 0, 0, s, t) 

+D{m\,m 2 Hi ,0,0,s,u)]j , (llf) 
where T| is the third component of the weak isospin of the quark running in the loop. 

III. NUMERICAL RESULTS 

We are now ready to illustrate the importance of squark loop contributions with a few 
examples. For the numerical analysis we have used the leading order CTEQ4L parameteriza- 



tion of the parton distribution function of the proton ||24|| , choosing the QCD renormalization 
and factorization scales to be the sum of the masses of the Higgs bosons in the final state. 
The effect of the running mass of the bottom quark can be very important, therefore, we 
have also included it in our calculations. In fact, when the the main contribution to the 
processes comes from bottom-quark loops and/or bottom-squark loops, the cross section 
is proportional to the <pbb Yukawa coupling to the fourth power. Taking 3 and 5 GeV for 
typical running and pole b— quarks mass respectively we can see that this effect can reduce 
the Higgs pair production by a factor (3/5) 4 ~ 1/8. 

As discussed in Sec. II. A. 2, the squark loop diagrams shown in Fig. 1 fall into three 
groups of diagrams, (3) + (4), (5) + (6), and (7) + (8), the sum of diagrams in each group 
being finite and gauge invariant. In unpolarized pp or pp scattering, where only the sum 
of the square of all helicity amplitudes is measurable, the squark contribution can therefore 
be characterized by three loop functions and the associated products of coupling constants. 
In order to assess the importance of these three sets of diagrams, we show in Fig. 3 their 
individual contributions to the subprocess cross section (a) for the production of hh pairs. 
For the sake of simplicity, only a single squark mass eigenstate (bi) was included here, 
whose mass is given on the x— axis. We chose = 100 GeV and = Sm^ = 300 GeV 
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as typical values. Note that diagrams (3), (4), (7) and (8) involve dimensionful couplings, 
while in diagrams (5) and (6) only dimensionless couplings appear. In order to show the 
corresponding loop functions, we have therefore considered fixed "typical" coupling constants 
(tan/5 = 50, Mg = 325 GeV, M A = 100 GeV, A t = A b = -410 GeV, [i =-640 GeV). 
However, these couplings were not varied as the mass of the squark in the loop is changed, 
while physical couplings do usually depend on the masses of the squark mass eigenstates, 
e.g. through the change of the q~L — qn mixing angle; see Appendix A. 

The loop function describing diagrams (3) + (4) is given in Eqs. (||a,b). It is the same 
(up to an overall factor) as that describing the squark loop contribution to single Higgs 
boson production ||. This contribution (dotted curve) seems to be important only when 
real H — > hh decays are possible and tan/? is not too large, since otherwise Br(H — > hh) 
becomes very small. The contribution from diagrams (5) + (6) (dashed curve), which involves 
quartic scalar couplings, is given in Eqs. (f|c,d). We find that this contribution can increase 
the total cross section by no more than a factor of a few. The reason is that this quartic 
scalar coupling cannot significantly exceed the square of the corresponding Yukawa coupling 
appearing in the quark loop contribution, and the squarks in the loop cannot be much lighter 
than the corresponding quarks. Of course, mg ^> m b , but < m t is still allowed. 

The potentially largest contribution therefore comes from diagrams (7) + (8) (solid line), 
Eqs. (|8]e-h), which involve trilinear Higgs-squark-squark couplings. These dimensionful 
couplings depend on unknown soft breaking parameters, and might be (much) larger than 
the mass of the lighter squark in the loop. On the other hand, while the loop functions 
for diagrams (3) + (4) and (5) + (6) slightly increase with increasing squark mass as long 
as rriq < \f§/2, the loop function for diagrams (7) + (8) starts to decrease as soon as 
rriq > m h /2 P5] . Once rriq > Vs/2, all squark loop functions become real, and drop rapidly 
with increasing rriq, approximately like m^ 4 . Due to the quick falling parton distribution 
functions, the largest contribution to the total Higgs pair production cross section come from 
values of s not far above threshold. Figure 3 then shows that squark loop contributions to 
the total cross section can only be large if the mass of the squarks in the loop does not much 
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exceed that of the produced Higgs bosons. 

The results of Fig. 3 allowed us to search the MSSM parameter space for parameters 
that maximize some subset of the diagrams listed in Fig. 1. In these searches we imposed 
the following constraints on parameters. First, most SUSY models predict 

1 < tan/3 < m ^ mt \ ~ 55 . (12) 
m b (m t ) 



Second, we have interpreted the unsuccessful search for Higgs bosons at LEP |26j to imply 
™>h > 90 GeV if the ZZh coupling has similar strength as the corresponding coupling in the 
SM. Otherwise the ZAh coupling is large, and (rriA + rrih) > 175 GeV is required. When 
computing the masses and couplings of the CP-even Higgs bosons, we have included squark 
and quark loop corrections |27]] as given by the 1-loop effective potential [p8|-|3~0"|l . Turning to 



the squark sector, for simplicity we took the same soft breaking mass rriq for m 4 ~ L = m~ bL , m- tR 
and ttii , and also took the same value A q for the soft breaking parameters A t and A b . The 
squark sector is then completely determined by rriq, A q , tan/3 and the supersymmetric Higgs 



mass parameter /i. In our scans we have imposed the LEP search limit |JT]] m^, m b > 80 
GeV *™ 



As anticipated, we found that the potentially largest squark loop contribution comes 
from diagrams (7) + (8). In Figs. 4a and 4b we show that (mostly) due to this contribution, 
squark loops can increase the total h— pair cross section at the Tevatron by more than 
two orders of magnitude if the mass of the lighter squark eigenstate in the loop is close 
to its experimental lower limit. In Fig. 4a we chose = 100 GeV and tan/3 = 50. The 
thick curves start at (m,, A q , = (325, —410, —640) GeV, corresponding to hb\bi coupling 
V h n = 455 GeV. This choice saturates the LEP Ah search limit (the ZZh coupling is very 
small here). The heavy solid and dashed curves have been obtained by varying rriq and 
fi, respectively, keeping all other parameters fixed. Note that increasing ii (i.e. decreasing 
\f/,\) not only increases m~ bi , but also reduces the hb\b\ coupling, leading to a rapid drop-off 
of the squark loop contribution. On the other hand, increasing rriq leads to a somewhat 
slower decrease of this coupling. In both cases the squark loop contribution to the total 
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cross section becomes essentially negligible for m~ bi > 150 GeV. For slightly smaller squark 
masses, there is mild destructive interference between quark and squark loops. 

In Fig. 4b we have chosen rriA = 500 GeV and tan/3 = 2.0. The thick curves originate 



at (rriq, A q , /j,) = (380,510, —975) GeV, which saturates the LEP ZZh search limit |33[] and 
gives V h n = 475 GeV. The thick heavy solid, dashed and dotted lines have been obtained 
by varying, one at a time, uiq, A q and /i respectively. We see that here the squark loop 
contribution remains significant out to ~ 200 GeV. This is partly due to the fact that 
the quark loop contribution in Fig. 4b is about fifty times smaller than in Fig. 4a, which in 
turn results from the large enhancement of the b— loop contribution compared to the SM, by 
roughly a factor tan 2 /3 = 2,500 in the amplitude. The bottom Yukawa coupling in Fig. 4a is 
nearly as large as the top Yukawa coupling in Fig. 4b. The former than gives a much larger 
quark loop contribution than the latter, since for the relevant values of s, the absolute value 
of the (mostly imaginary) b— quark loop function is much larger than that of the (mostly 
real) t— quark loop function. 

This enhancement of the contribution of b— quark loops at large tan/3 also implies that 
prospects for detecting a signal for the production of neutral Higgs boson pairs at the next 
Tevatron collider run might be better if tan/3 is large. We estimate that a cross section of 50 
fb or more might be detectable. This would lead to roughly 10 events per experiment, each 
with 4 high— tagged 6— jets and double bb invariant mass peaking, assuming an integrated 
luminosity of 2 fb -1 and an overall efficiency of 10%. The enhancement of the pure quark 
loop contribution required to achieve this cross section is given by the thin lines in Figs. 4a 
and 4b. We see that for tan/3 = 50 the total cross section can exceed this sensitivity limit 
by more than an order of magnitude. However, the maximal light squark mass compatible 
with such a large hh production cross section is about the same at low and at large tan/? 
(~ 110 GeV). 

It is interesting to notice that in Fig. 4a, rrih is significantly below tua, i.e. the difference 
is up to 25 GeV, for all values of mg where squark loop contributions are significant. This 
mass pattern, which is quite unusual for large tan/?, is due to non-logarithmic corrections 
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to the Higgs boson mass matrix involving trilinear scalar interactions; technically, due to 



large contributions to A 12 , in the notation of refs. [£8^9]. This means that the total cross 
sections for the production of pairs of other Higgs bosons are quite small at the Tevatron. 

At high tan/3 the Higgs pair production can be enhanced by factors tan 4 (3 or 1/ cos 4 (3 
with respect to the SM production mechanism. The only exception to this is the hh {HH) 
channel in the large (small) Ma limit, in which the factor sin a/ cos (3 (cos a/ cos (3) appearing 
in Yukawa couplings goes to 1. For instance, in the case that (nig, A q , ji) = (1, 1, 1) TeV, 
tan/3 = 50, and M^=100 GeV the total cross section at the Tevatron for the hh, HH, and 
AA production is 3.3, 0.034, 3.9 fb respectively. When Ma is increased to 130 GeV these 
cross sections change to 0.15, 0.1, and 0.7 fb. 

The Higgs pair production cross sections do remain sizable at the LHC In Figs. 5a-d 
we show the squark loop contribution to hh (a), HH (b), hH (c), and AA (d) productions. 
The parameters taken in Figs. 5a-c are the same as in Fig. 4a, which had been chosen to 
maximize the hbibi coupling. The comparison of Figs. 4a and 5a shows that the relative 
importance of squark loops is almost independent of the center-of-mass energy. Of course, 
the total cross section increases greatly when going from y/s = 2 TeV to 14 TeV due to 
the rapid increase in the gluon-gluon luminosity. This is illustrated by the thin lines, which 
again correspond to a total cross section of 50 fb. The quark (mostly b) loop contribution 
by itself now exceeds this cross section. Nevertheless, the background also increases when 
going from the Tevatron to the LHC. In the absence of a dedicated analysis of signal and 
background, we do not want to claim that a total cross section of 50 fb necessarily gives a 
detectable signal at the LHC, in spite of its considerably higher anticipated luminosity. 

Although the starting point of the curves in Figs. 5a-c had been chosen to maximize 
^/ibibi' we fi n d a very large squark loop contribution also for HH production (Fig. 5b). In 
this case the squark loop contribution at first increases with increasing rag . The reason 
is that the Hb\b\ coupling increases quickly, from ~ 150 GeV at the starting point of the 
curves to ~ 350 GeV near the maximum of the dark solid curve. In fact, this coupling keeps 
increasing even further as rriq is increased, eventually reaching ~ 450 GeV. However, this 
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increase is overpowered by the rapid drop of the loop function once m~ bi significantly exceeds 
m H (see Fig. 3). 

The biggest relative contribution from squark loops appears in hH production, Fig. 5c, 
giving rise to an enhancement factor ~ 500 in some cases. In this region of parameter 
space the Hbb Yukawa coupling becomes very small, due to the unusual mixing pattern of 
CP-even Higgs bosons caused by radiative corrections at large tan/3 and large \A q \ and 
At the peak of the curves the Hbb coupling vanishes completely, and the hti coupling is 
quite small, leading to a very small total quark loop contribution (see the behavior of the 
thin lines). Notwithstanding, when m q or jj, are increased beyond this point, the quark loop 
contribution reasserts itself while the squark loop contribution decreases in absolute size, 
leading to a steep drop of the relative importance of the squark contribution: it becomes 
essentially negligible for > 220 GeV. The suppression of the Hbb coupling also explains 
why squark loop contributions to HH production can remain significant up to m~ bi ~ 300 
GeV (see Fig. 5b). In this case the dominant quark contribution comes from top quark 
loops, so squark loops remain significant unless m~ raf . 

For the parameters chosen in Figs. 5a-c the squark loop contribution to AA production 
is totally negligible. In this case diagrams (7) and (8) in Fig. 1 have to include at least 
one heavy squark mass eigenstate, since the diagonal Abib\ and At\t\ couplings vanish 
identically. The off-diagonal AbJ>2 coupling is actually quite large, ~ 270 GeV at the 
starting point of the curves in Fig. 5a. However, at the same time m~ h2 = 455 GeV, which 
suppresses the contributions from diagrams (7) + (8) to an insignificant level. 

In Fig. 5d, we therefore show results for a scenario with relatively small my. rriA = 150 
GeV, tan/5 = 4, A q = -110 GeV, p: = 345 GeV, and m q between 115 and 350 GeV. The same 
set of parameters also yields a relatively large quartic AAtit\ coupling, so that diagrams 
(5) + (6) in Fig. 1 are maximized. Since tan/3 is fairly small, the quark loop contribution is 
dominated by top quark loops; however, the loop function is also different for CP-odd Higgs 
bosons, see Eqs.(6), and its contribution is suppressed by a factor cot 4 /? ~ 1/250 compared 
to the SM case. Therefore diagrams (7) + (8) can make a large relative contribution, as 
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long as rriq < m t . Even though squark loop contributions can increase the AA cross section 
by a factor of about 200 in this case, it still remains well below the cross section for HH 
production. In scenarios where the AA cross section is comparable to the HH cross section, 
we find that squark loop corrections to AA production are quite modest. 

Given that we are only working in leading order in QCD, our predictions for the total 
cross sections have significant uncertainties due to the choice of scale in a s , my,, and the 
parton distribution functions. We took the same scale everywhere, viz. the sum of the 
masses of the produced Higgs bosons. One could therefore infer the existence of squark loop 
contributions to the total cross section only if it changes the quark loop result by at least 



a factor of two ||34j| . However, smaller squark loop contributions might still be visible in 
some distributions. For example, for large tan/3 the quark loop contribution is dominated 
by b— quark loops (except for hh production at large m^, where h is always SM-like). Since 
mg 3> m&, the squark and quark loop functions show quite a different dependence on s. 
This is illustrated in Fig. 6, where we show da(hh)/dy§ for the point in Fig. 5a where b 
and b loops contribute equally. The solid histogram shows the contribution of b loops only, 
while the dashed histogram includes squark loops; the total cross section (area under the 
curves) has been normalized to be the same in both cases. The b loop contribution peaks 
just beyond the threshold at \/§ = 2m/,, but the b\ loop contribution clearly shows up as 
a second peak at vl ~ 2m^. This distribution is in principle directly measurable if both 
Higgs bosons decay hadronically (with combined branching ratio ~ 80%). Given sufficient 
statistics, one might be able to infer the squark loop contribution in this manner even if the 
total cross section is dominated by quark loops. 



IV. SUMMARY AND CONCLUSIONS 

In this paper we have calculated squark loop contributions to the pair production of 
two neutral Higgs bosons. If CP is conserved, squark loops contribute only if the two 
produced Higgs bosons have identical CP quantum numbers. In Sec. II we gave complete 
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analytical expressions that allow the evaluation of these contributions. For completeness we 



also included expressions for the quark loop contributions |?0] . These formulae are written in 
a completely general fashion with the explicit expressions for the relevant coupling constants 
in the framework of the MSSM being listed in Appendix A. 

In Sec. Ill we showed some numerical results for the MSSM. We found that squark 
loop contributions can increase the total cross section for the production of two CP-even 
Higgs bosons by more than two orders of magnitude. However, such large contributions are 
possible only if three conditions are satisfied: 

• The relevant Yukawa coupling must be large. In case of the (s)top, this is always 
true for at least one of the two CP-even Higgs bosons of the MSSM. However, as 



pointed out quite some time ago [35], in the MSSM (as in other models with more 
than one Higgs doublet) the bottom Yukawa coupling can also be large, if the vacuum 
expectation value that gives rise to the mass of the b quark is small. In the MSSM 
this happens for large values of tan/?. However, the bottom Yukawa coupling is not 
expected to exceed that of the top quark. 

The lighter of the two superpartner whose Yukawa coupling is large (generally b\ at 
large tan/? and t\ at small tan/3) must not be much heavier than the Higgs bosons in 
the final state. This condition is especially critical at large tan/3, since here b\ loops 
have to compete with b quark loops. For equal Yukawa couplings (i.e. tan/3 ~ m t /m^j, 
the squared b— loop contribution exceeds the squared t— loop contribution by a factor 
of ~ 50. If this condition is satisfied, open squark pair production should be detectable 
at the same collider, unless the squark-LSP mass difference is very small. 

The relevant trilinear soft breaking parameters and/or must be significantly larger 
than the mass of the lighter squark eigenstate. This implies that the diagonal entries 
of the corresponding squark mass matrix must also significantly exceed the smaller 
eigenvalue of this matrix. This requires a modest amount of fine tuning. However, in 



18 



the absence of a complete theory of Supersymmetry breaking this possibility should 
not be discounted. 



If these three conditions are satisfied, the region uia — 100 GeV can perhaps even 
be probed at the next run of the Tevatron collider, using searches for final states with 
4 high— px b— jets. Moreover, this process has a big advantage over the <pbb associated 
production |TJ|, which leads to the same final state: the reconstruction efficiency of Higgs 
pairs is expected to be one order higher than that for <pbb production. This originates from 
the fact that all 6-jets coming from Higgs boson pair decays are energetic, while the two 
associated 6-jets in the <frbb channel are quite soft. At the same time the process under study 
has a cross section not much smaller than that for (pbb. The reach of the LHC should be 
much higher, but a quantitative statement will only be possible after a detailed analysis of 
signal and background |36| . Furthermore, we should point out that our result for the 1-loop 
cross sections are probably very conservative, since we can expect large QCD corrections, 
which could increase the cross section by as much as a factor of ~ 2 |7|,|§|,[2~1"|1 . 

Squark loop contributions to the pair production of two CP-odd Higgs bosons are more 
modest in general. In this case trilinear Aq^q* couplings contribute only if at least one 
squark in the loop is a heavy mass eigenstate, which leads to a suppression of the squark 
contributions. Nevertheless, squark loops can give rise to large enhancements of the A pair 
production cross section if < m t and tan/3 ~ 5, which leads to a very small quark loop 
contribution. However, the total cross section for AA production remains quite small in this 
case. 

In general AA final states, as well as the HA and HA final states, which receive no 
contributions from squark loops (but do receive Drell-Yan like contributions from light qq 



annihilation [20]), can be significant. Note that in most regions of parameter space the CP- 
odd Higgs boson A is nearly degenerate with one of the two CP-even Higgs bosons. This 
means that often three different channels {e.g. HH, HA and AA) contribute to essentially 
the same final state, and therefore have to be added. This obviously increases the chance to 
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detect a signal for Higgs pair production in the 46 channel. At the same time it complicates 
the interpretation of such a signal, e.g. the extraction of the relevant coupling constants. 
Nevertheless, we saw in Fig. 6 that the analysis of the 46 invariant mass distribution can 
help to disentangle the various contributions to the signal. We are therefore hopeful that 
the search for the pair production of neutral Higgs bosons in the 46 channel will provide 
information that will help us to pin down the Higgs sector, and perhaps also the squark 
sector of the theory. 
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APPENDIX A: MSSM COUPLING CONSTANTS 



We denote the weak mixing angle and couplings as 

s w = sm9 w , c w = cos9 w , g = e/s w , g z = g/c w 



We define squark mass eigenstates via 



7r 




cos9f sin6*/ ' 




7l' 


. f'2 . 




k — sin9f cosOf . 







(Al) 



wnere uie mixing angle uj is aennea in sue 

with eigenvalues m 2 ? and m 2 ? (m? < mi ): 

to /1 h v n H> 



, f = u,d . (A2) 
t the mass matrix becomes diagonal 
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with 



m ) L =m h + m ) + M l cos 2 P( T i - Q f s w) ; ( A4a ) 

m 2 j R = mj R + mj + M§ cos 2/3Q /s^ ; (A4b) 
m u (ii cot P - A* u ), f = u 

(A4c) 

ra/(/itan/3 - A* f ), f = d 



ml = < 



Squark current eigenstate bilinears can then be expressed in terms of current mass as 
follows: 

flh = +c)flh - CfSffif 2 - Cfsjfh + s)fff 2 - (A5a) 
Mr = +cfffifi + cj/r/a - - CfSffif*, (A5b) 

/£/l = +c fSf h~h - s)fth + - c /S/ /|/ 2 ; (A5c) 

/S/r = +4/1/1 + C /S/ /f / 2 + C/S//I/! + 4/|/ 2 , (A5d) 

with Cf = cos9j and Sf — sinfy- for f = u,d. 

After these preliminaries we are ready to list the relevant couplings. We list couplings to 
squark current eigenstates only; these can be converted into couplings to mass eigenstates 
using Eqs. (17). 

H(h, A)-Fermion-Fermion 

uuh° : - Z uC ° Sa « uuH° : - 9 g ™" sin " uuA° : cot ^75 

2Mw sin p 2Mw sin p 2Mw " ' ° 

ddh° : + 9 y sin ° : - Z dC ° Sa R ddA° : +^^tan/5 75 

2Mwcosp 2Mwcosp 2Mw 

Gluon— Squark— Squark 

where k and p are the incoming momenta of the squarks g«, with % = 1,2. 
H-Squark-Squark 
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u\ u L H° : 
u* r urH° : 
u* L u R H°,u* R u L H° : 
d* L d L H° : 
d* R d R H° : 
d* L d R H°,d R d L H : 
h-Squark-Squark 

hf*f = Hf*f{sina — > cos a, cos a — > — sina,sin(a + /3) — > cos(a + (3), 
cos(a + /?) — > — sin(a + /?)} 

A-Squark-Squark 

u* L u R A°: -iJ^(A u cot(3-ii) u* R u L A° 

d* L d R A°: -i^AatwP-v) d* R d L A° 
Gluon-Gluon-Squark-Squark 

y^y yj y« • y s 2 2 



H-H-Squark-Squark 



^ L (iJ ) 2 : 


[g 2 »Tt 2 sin 2 a 
UM^sin 2 /3 


+ \g\{Tzu - 




cos 2a] 


u R u R (Hy : 


f <? ?n 2 sin a 




cos 2a] 




UM^sin 2 /3 




dld L (H ) 2 : 


|-g 2 m^ cos 2 a 
Um^cos 2 /3 


+ \g 2 z {T, d - 


- sfvQd) 


cos 2a] 


d* R d R (H ) 2 : 


r <y 2 m 2 cos 2 a 


+ y z s wQd 


cos 2a] 




Um^cos 2 /? 





h-h-Squark-Squark 



hhqq = HHqq{sm a — > — cos a, cos a — > sin a, cos 2a — > — cos 2a} 
H-h-Squark-Squark 

Hhqq = HHqq{cos 2a — > — sin 2a, sin 2 a — > sin 2a, cos 2 a — > — sin 2a} 
A-A-Squark-Squark 



-(^ + l9zMzcos(a + P)s^ 

^ 2M^np ( Au sma + /icosa) 
-(^S? + ^M z oos(a + + 
-( 9 ^-y Z M z cos(a + ^ w ) 

^2l^^(^ COSa + / iSina ) 



+i^(A u cot/3-//) 
+ i^-(A d tan/3-/,) 
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AAqq = HHqq{sin a — > — cos (3, cos a — > sin /3, cos 2a — > — cos 2a} 

Higgs°-Higgs°-Higgs° 

(H ) 3 : -2fM z cos2acos(a + /3) (/i ) 3 : -2fM z cos 2a sin(a + /3) 
(A°) 2 H° : +2fM z cos2 / 3cos(a + / 3) (A°) 2 /i° : - 9 -f M z cos2/3sin(a + /3) 
(/i°) 2 #° : -^M z [2sin2asin(a + /3) -cos2acos(a + /3)] 
(#°) 2 /i° : + 2fM z [2sin2acos(a + / 3) + cos2asin(a + /?)] 
Z-Higgs-Higgs 

h°{p)A°{k)Z fl : +is*cos(a-/?)(fc-p)" /i (p)A (A;)Z M : sin(a - /3)(Jfe - p)» 

APPENDIX B: CONVENTIONS 

We chose the following polarization vectors in our calculations 

4 = <^2 = ^=(o,-i,-«,o) , 
e r = 4 = ^(M,-*,o) , 

where the first gluon is moving along the z axis. 

We have used scalar C and D Passarino-Veltman functions in our analytical formulas 
which are expressed as: 

C \plO,pl2, P 20,mO,ml,m2] = 

(m 2 )- 1 j d'q([q 2 - ml] [(q + Pl f - m\] [(q + p 2 f - m,])' 1 , (Bla) 
D [pl0,pl2,p23,p30,p20,pl3, mO, ml, m2, m3] = 

(tTT 2 )- 1 J d 4 q([q 2 - m 2 ] [(q + Pl f - m\] [(q + p 2 f - m 2 ] [(q + p 3 f - m 2 ])' 1 . (Bib) 
Our convention for the scalar arguments is piO = pf, pij = (p, — Pj) 2 and mi = m^. 
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FIG. 1. Feynman diagrams for hh, HH, hH, and AA Higgs boson pair production. = h,H 

for = 1,2 respectively, = 91,92 for = 1,2. The crossed diagrams are not shown. 
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FIG. 2. Feynman diagrams for the hA and HA Higgs boson pair production. Hj = h,H 
1,2 respectively, q^m = <7i,<?2 for k(l) = 1,2. The crossed diagrams are not shown. 
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FIG. 3. Contributions of the (3) + (4), (5) + (6), and (7) + (8) diagram sets to the subprocess 
cross section for the production of hh pairs as a function of . We included only the b± effects 
and assumed nth = 100 GeV and \/i = 300 GeV. 
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FIG. 4. Ratio of the total cross section to the cross section including only the quark contribu- 
tions for gg — > hh at the Tevatron. (a) corresponds to a large tan/3 (= 50) scenario and we took 
(rriq, A q , fi) = (325, —410, —640) GeV. The heavy solid and dashed curves have been obtained by 
varying, one at a time, rriq and fj,, respectively. In (b) we display a low tan/3 (= 2) scenario using 
(mq,A q ,fi) = (380,510,-975) GeV and conventions as in (a). Here, the heavy dotted line was 
obtained by varying A q . In both (a) and (b), the thin lines correspond to the enhancement needed 
for the total cross section to be at the level of 50 fb. 
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FIG. 5. Same as in Fig. | for the hh (a), HH (b), Hh (c) and AA (d) production at the LHC. 
The conventions are the same as in Fig. In (a)— (c) we chose the parameters used for Fig. 4a. In 
(d) we assumed that tua = 150 GeV, tan/? = 4, A q = — 110 GeV, fi = 345 GeV, and niq between 
115 and 350 GeV. 
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FIG. 6. do jd\Tk including only quark loops (solid line) and considering quark and squark loops 
(dashed line). We have chosen (m ? -,Aj,/i) = (335,-410 - 640), M A = 100 GeV, and tan/3 = 50 
for which quark and squark loops contribute almost equally: Ototal = 10. fb, a squa rk = 5.6 fb. 
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